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ABSTRACT 

We use Smoothed Particle Hydrodynamics to simulate the formation of a massive 
(IO^Mq) stellar cluster system formed from the gravitational collapse of a turbulent 
molecular cloud. We investigate the hierarchical clustering properties of our model 
system and we study the influence of the photoionizing radiation produced by the 
system's multiple 0-type stars on the evolution of the protocluster. We find that 
dense gas near the ionizing sources prevents the radiation from eroding the filaments 
in which most of the star formation occurs and that instead, ionized gas fills pre- 
existing voids and bubbles originally created by the turbulent velocity field. 
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1 INTRODUCTION 

Massive clusters and OB associations associat ed with giant 
HII r e gions such as R136 iii 30 Doradus (e.g. Selman et al.l 
I1999I : iTownslev et"all I2OO6I : llndebetouw et all |2009|) the 
central clusters i n NGC 3603 (e.g. iTapia et aLfbOOll: 
Stolte et all l2006h and NGC 604 (e.g iTiillmann et all 



dynamical effect on clusters of removing gas left o ver b y 



the star f ormation process, while, e.g, 
iTenorio-Tagle fc Bodenheimeil jigS^ I. [Franco et atl (|l994l 'l 



LarsonI (Il987 



20081). and Trump l er 14 15 and 16 and Collinder 228 (e.g. 
Tapia et all I2OO3I '2006") in the Carina HII region, each 
containing tens to hundreds of 0-type stars, are ideal 
places to test our understanding of star formation and the 
effects of feedback from the high-mass end of the IMF. 
Star formation in such regions appears to be hierarchical, 
with structure extant over a wide range of length and 
mass scales. The dynamics of these regions, or at least 
of their remaining gas content, seem to be dominated by 
feedback from O-stars in the form of giant or multiple HII 
regions, wind-blown bubbles and supernova shells. These 
phenomena drive powerful outflows, expelling mass from the 
embedded clu sters and perhaps unbinding them. Several 
authors (e.g iLada fc Ladal l2003l : iFall et all l2005l ) have 
concluded that the majority of star clusters do not survive 
their embedded phase and the loss of large quantities of 
gas via outflows offers an attractive explanation for this 
phenomenon. 

The dispersal of clouds by t he ac t ion o f feedback 
has been studie d by many auth ors. Hills' (l9 80|) and later 
iGoodwinI ll 19971). iBoilv fc Kroupa (2003b). B oilv fc Kroupal 
l|2003al ) and iGoodwin fc Bastieml (fiood ) considered the 
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and iMatzneJ (|2002l ) consider the self-limiting effect of 
feedback on the star-formation efficiency in clouds. These 
authors did not, however, mo d el the process of gas re- 
moval itself. iGever fc Burker^ (|200ll ) performed hybrid 
N-body/SPH calculations in which they modelled the 
interaction of the gas ad stars self-consistently, but the 
expulsion of the gas, by heating or by driving a Shockwave 
through it, were still parameterized. It is clear that an 
understanding of the process by which stellar systems lose 

gas is crucial. 

iTenorio-Tagld l|l979l ) introduced the concept to the 
champagne flow, in which an O-star born near the edge of 
a cloud creates a cavity around itself through which gas 
eroded from the wall of the cloud flows at high veloc ity into 
the ISM. Using an analytic model. Iwhitworthi (|l979l ) found 
that this process would result in star formation efflciencies 
at the scale of whole giant molec ular clo uds wou ld be ~_5%. 
2-D numerical simulations by iBodenh cimer et al.l (|l979l ) 
also showed that champagne flows could be an efficient 
dispersal mechanism and they c learly do play such a role 
in m any systems, such as O rion ([Hennev et al. l2005h . How- 
ever, iMazurek* (1980) and lYorke et al.l l|l989t ) found that 
the disruption of clouds by HII regions is strongly impeded 
or prevented if the clouds are in freefall and if the ionizing 
sources are located near the centres of the clouds. As noted 
above, observations often indicate that massive stars are 
more likely to be found near the centres of star-forming 
regions, or at any rate at the centres of subclusters, and 
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many simulations reproduce this r esult. This can be the 
result of competitive accre tion (e.R Bonnell et al.ll200^ ) or 
rapid mass segregation (e.g lMcMillan et al.ll200'^ 

Additionally, the above-cited models of champagne 
flows treat the gas as being initially uniform or at least 
smoo thly-varying. In contrast, simulations and observations 
(e.g. iLongmore et al. 

I I2OO7) reveal that massive stars are 
often embedded in highly inhomogeneous gas from which 
they are continually accreting and which can se verely 
limit the effects of feedback (e.g IPeters et al.ll2010bl '). The 
com plex stru cture of the gas may derive from turbulence 
(e.g iKlessenI [2OOII. or from gravitational amplification of 
even very small quantities of noise in the density field 
l|BonneU fc Batd [iooaV The effects of feedback on such 
systems demands the use of 3-D sel f-gravitating simu la- 
tions such as those c arried out by, e.g, iDale et al. 1 (120051 ) or 
IPeters et all (|2010bh . lOale et al.1 (|2005h found that strongly 
inhomogenous and anisotropic accretion fiows onto an 
ionizing source strongly limit the effect of ionizing feedback, 
with radiation only being able to escape the immediate 
vicinity of the source over a limited fraction of the sky, in 
strong contrast to 1-D models of HII regions as Stromgren 
spheres. 

Feedback is also thought to be responsible for trig- 
gering additional star formation in m any regions an d 



there is much literature on this topic (e.g Dcharve ng et al 



20031 : IZavagno et all I2OO6I : iKoenig et all |2008; P uga et al 



20091 ). However, identifying triggered star formation from 
an observational perspective is extremely difficult and 
morphological indicators, such as shell-like structures of 
stars and gas, or the close proximity of young stellar objects 
to HII regions must often be relied upon. 

In this paper, we investigate the evolution of a massive 
cluster (or, more accurately, a massive complex of clusters) 
formed from the gravitational fragmentation of a lO^'M© 
turbulent molecular cloud. We examine the clustering 
properties of the stellar systems formed in the cloud and we 
include the effects of the photoionizing radiation from the 
multiple O-type stars hosted by the clusters on the highly 
inhomogeneous gas remaining in the system as a first step in 
understanding the transition between large-scale embedded 
star formation and open clusters or OB associations. Our 
numerical methods are presented in Section 2. Section 3 
contains our results, which are then discussed in Section 4. 
Section 5 contains our conclusions. 



2 NUMERICAL METHODS 

We use the w ell-known Smoot hed Particle Hydrodynam- 
ics technique (|Monaghanl Il992l ) , which is ideal for study- 
ing the evolution of molecular clouds and embedded clus- 
ters. We use a version of the Benz code (|Benzll 19901 ) and 
the standard artificial viscosity prescription, with a = 1, 
P = 2. Particles are evolved on individual timesteps. The 
code is a hybrid N-body SPH code in w hich stars are rep - 
resented by point-mass sink particles (|Bate et al.l Il995h . 
Self-gravitational forces between gas particles are calculated 
using a binary tree, whereas gravitational forces involving 
sink-particles are computed by direct summation. Sink par- 
ticles are formed dynamically and may accrete gas particles 



and grow in mass. 

The thermodynamics of the gas are treated using a 
modified Larson-type piecewise equation of state (|Larsonl 
i2005,) comprised of three barotropic equations of state 
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tion of state mimics the effects of line cooling. The 7 = 1.0 
segment approximates the effect of dust cooling whilst 
the 7 = 1.4 segment mimics the effects of when collapsing 
cores become optically thick to infrared radiation. The final 
isothermal phase of the equation of state is simply in order 
to allow sink-particle formation to occur, which requires a 
subvirial collapsing fragment. 

In these simulations, the model cloud has a mass of 
IO^Mq and is seeded with a Kolmogorov turbulent velocity 
field whose total kinetic energy is equal in magnitude to 
the cloud's initial gravitational binding energy, so that the 
cloud is marginally bound and is therefore expected to form 
stars rapidly and efficiently in the absence of feedback. The 
initial mean Jeans mass is ~ 260 M0, so that the initial 
conditions contain ~ 4000 thermal Jeans masses and we 
would therefore expect of order this number of fragments 
to form if the system were isothermal. However, the mean 
intitial temperature is ~ 40K, and the collapsing regions 
cool to lOK, and the Jeans mass in the dense filaments that 
form is therefore lower. 

We model the cloud at two resolutions, using 2.5x10^ 
and 10^ particles. We use the higher resolution calculation 
to examine the clustering properties of the star formation 
in the system, and the low-resolution run to study the 
effects of ionization. Due to the large number of ionizing 
sources, we found that simulating ionization in the high- 
resolution run was too time-consuming. M ass resolution in 
the h igh-resolution calculation was 4M0 (|Bate fc BurkertI 
Il997h . but in the low-resolution simulation is rather large 
at IOOMq. The sink particles formed in this simulation 
cannot therefore be regarded as stars, but are instead 
small clusters. If such a cluster is sufficiently massive, it 
will host one or more 0-stars. For all sink particles whose 
mass exceeds 600 Mq we compute the number N30 of 30 
M0 O-stars assuming that the clusters have a Salpeter 
IMF with minimum and maximum stellar masses of O.IMq 
and IOOMq. We then take the ionizing luminosity of each 
cluster to be N30 x 10*^ ionizing photons per second. 

The effects of ionization on the gas are simulated using 
a modified version of the S tromgren volume algorithm 
described in lDale et al.l (|2007l ) to identify ionized particles, 
whose temperatures are then set to lO^K and treated 
isothermally. The algorithm has been modified in two 
ways. Firstly, the original algorithm updates the ionization 
fraction of all particles on the shortest dynamical timestep 
in use by any particle in the code. Instead, the new code 
only updates the ionization states of particles when they are 
being dynamically updated, so that the ionization timestep 
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of each particle is taken to be equal to its dynamical 
timestep. 

Secondly, it was necessary to modify the ionization 
algorithm to allow simulation of the effects of multiple 
ionizing sources. In the case that no two HII regions overlap, 
no changes to the algorithm would be required, but if the 
ionized regions excited by two or more different sources 
overlap, a given particle may be receiving flux from more 
than one source. For an ionized particle, the recombination 
rate per unit volume of the ionized gas making up the 
particle is aBU^, where n is the number density of ions 
and Ob is the case-B recombination coefficient. To keep 
the particle ionized, ionizing photons must be subtracted 
from radiation beams passing through the particle. The 
problem that must be solved is to determine how many 
photons should be subtracted from each beam. We solve 
this problem in a crude way by assuming that all sources 
illuminating a given ionized particle are equally responsible 
for keeping it ionized. Strictly, the global ionization state of 
the whole cloud should then be found by iteration, but this 
is time-consuming. Instead, at a given timestep, if a gas 
particle is receiving radiation from A'^ ionizing sources, when 
the particle's ionization state is next updated, we assume 
that the volume recombination rate subtracted from the 
radiation beam of each source passing through the particle 
is oiBn^ /N. We therefore assume that the recombination 
load is shared equally among all the sources illuminating a 
given particle. 



3 RESULTS 
3.1 Clustering 

Both the low resolution and the high resolution clouds evolve 
in a similar fashion over the 4.3 Myr that the evolutions 
were followed (without feedback) . The turbulence in the gas 
generates significant internal structure in the form of fil- 
aments. Regions within the filamentary structure become 
self-gravitating and collapse to be replaced by sink parti- 
cles (which represent stars in the high-resolution simulation, 
and small clusters in the low-resolution run). Larger-scale 
regions become gravitationally bound due to the shock dissi- 
pation of kinetic energy, allowing clusters or groups of clus- 
ters to form at the junctions of the filaments. These forming 
clusters continue to grow by accreting infalling gas and other 
sink particles. The clustering properties of the model are of 
intrinsic interest for comparison to real massive star form- 
ing regions, but also since the distribution of clusters and 
massive stars determines to a large extent the infiuence of 
feedback. 

The low resolution run created 601 sink-particles 
whereas the high-resolution simulation created 5750 sinks. 
In Figure [1] we plot the location of all stars at a time of 
4.3Myr in this simulation. Amongst the sinks are several 
very significant clusterings including the most massive clus- 
ter with a mass density in excess of lO'^ M0 pc-3 inside 0.25 
pc. Figures [2] and [3] show the mass density for each sink in 
the high-resolution computation computed from the mass 
within Ipc of that sink, and from the distance to its tenth 
nearest neighbouring sink respectively. We see from these 
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Figure 1. Locations of sink particles in the high— resolution cal- 
culation. Coloured symbols indicate the centres of clusters whose 
densities exceed 10^ (blue), 10^ (green) and 10^ (red) Mq pc~^. 
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Figure 2. Plot of cluster density (as computed by the mass within 
Ipc) against cluster mass. 



figures that most of the more massive sinks are in highly 
clustered regions with local mass densities reaching 10* 
pc~^. In addition to the high clustering of most of the mas- 
sive sinks, we see a number of relatively massive sinks with 
masses ^ 100 M0 that are in low density regions as mea- 
sured from the distance to their tenth nearest neighbour. 
Such sinks likely represent small clusters that may contain 
what appear to be isolated massive stars in the context of 
the whole system. 



4 



James E. Dale, Ian Bonnell 



10" 
lO'' 
10* 
10" 
10* 
1000 
100 
10 

1 

0,1 
0.01 
0.001 




100 1000 10* 

Mass density (in 0.25pc) 



10° 



a 




Figure 3. Plot of cluster density (as computed by finding the 
volume containing 10 stars) against cluster mass. 



Figure 5. Positions of clusters (those massive enough to be ion- 
izing sources are marked as crosses) at the point when ionization 
was turned on. 




Figure 4. Column-density map viewed down the z-axis of the 
system at the point where ionizing radiation was turned on. 



3.2 Effects of ionization feedback 

We perform a separate low-resolution calculation in which 
we enable ionization when the system has an age of 3.7 Myr 
(0.64 tfi), which is approximately the time when the system 
begins to give rise to clusters large enough to host O-stars. 
At this point in the system's evolution, there are eleven such 
clusters and these therefore form our ionizing sources. Figure 
[4] shows a column-density plot as viewed down the z-axis 
of the system at this epoch, with white dots representing 



clusters. The gas distribution is composed of a network of 
dense filaments in which most of the star formation is occur- 
ring, particularly at their junctions, permeated by a froth of 
low-density voids. In Figure [5] plotted to the same linear 
dimensions, we show the positions of the ionizing sources 
(obviously as time progresses and clusters accrete, more ob- 
jects acquire enough mass to become ionizing sources). 

In Figure [6l we show at three different epochs column 
density maps of the cold gas (red in the left and centre panels 
) and the HII (blue in the right and centre panels). It is clear 
that the hot, low-density HII fills up cavities in the cold gas, 
so that the HII bubbles are illuminated from the edge rather 
than the centre, often by several widely-separated O-stars, 
and that several small HII regions eventually join to form a 
single large and irregularly-shaped one. The brightest part 
of the HII region is not centred on the main concentration 
of clusters, but instead abuts it. Additionally, as shown in 
Figure [J] these cavities existed before the ionizing sources 
were turned on and were generated by the turbulent veloc- 
ity field with which the cloud was seeded, so they have not 
been created by the action of feedback. We also note that 
the boundaries of the cavities are delineated by the dense 
filaments in which star formation is taking place, so that 
the morphology of the system appears to consist of several 
overlapping Hll-filled bubbles with vigorous s tar formation 
activ i ty on their border s. Several authors (e.g. lKoenig et all 
l2008l : iPuga et a l."200Sf) have taken this configuration to be 
indicative of triggered star formation but in our calculations, 
this is clearly not so. The overall HII morphology is irreg- 
ular, with an approximately round core centred at about (- 
10,-10) pc and approximately 25 pc across, and several spurs 
and offshoots, one extending about 50 pc towards the top 
right in the final frame of the simulation. The distribution 
of HII is also markedly clumpy due to the dumpiness of the 
gas and the action of multiple ionizing sources. In the final 
frame of Figure |6l it is clear that radiation and ionized gas 
are leaking out of the cluster through holes in the neutral 
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material, particularly in the bottom left of Figure 6c. This 
structure superficially resembles a champagne flow, but its 
origin is clearly rather different, since the bubble was not 
created by the HII region. 

We find, as in iDale et all l|2005l ). that the effect of the 
ionization is strongly constrained by the dense filamentary 
gas in which all the ionizing sources are embedded, and on 
which it has only a small effect. Ionization fronts propagate 
only into the low-density gas, so that only a small amount 
of dense material is ionized, but most of the HII is gener- 
ated from lower-density gas. In Figure [9l we plot the final 
density versus the initial density for all ionized gas (red) 
and for a randomly-chosen one percent of the neutral gas. 
In general, the density of the neutral gas has increased dur- 
ing the simulation, since the cloud is collapsing and forming 
stars. Conversely, the density of ionized gas has generally 
decreased (although not in all cases, as some material has 
become denser than it was initially before being ionized). 
In Figure IIOI we plot the cumulative mass of ionized gas 
against its initial density. The ionized gas is clearly drawn 
from material that was of a low density initially and much 
of it has subsequently become less dense as it expands into 
the voids between the filaments. 

The reason that photoionization has such a limited ef- 
fect on the dense gas in the cluster is that all the ionizing 
sources are embedded in, and continuously accreting from, 
dense material and most of their photon fluxes are absorbed 
by the accretion flows, leading to evolution ver y different 
from the classical picture of a champagne flow. IWalmslevI 
l|l995l ) gave a simple formula for a critical accretion rate 
Merit onto an ionizing source of mass M and ionizing lu- 
minosity L (Lyman continuum photons per second) above 
which all the ionizing photons are absorbed by the accretion 
flow, so that the HII region cannot expand and may instead 
contract: 

^->(^^^^)*. (3) 

We find that for all our sources, at virtually all times, 
their accretion rates are supercritical. If the accretion flows 
were spherically symmetric, they would then absorb all 
of the combined photon flux of the sources in the cluster. 
There are three ways in which such an accretion-swamped 
ionizing source may still grow an HII region. The first 
and most obvious is that the accretion rate may decline. 
Although the accretion flows in our simulation are by no 
means steady, the accretion rates do not decrease with 
time. The second possibility is that, as the source accretes 
and grows in mass, its ionizing luminosity increases as a 
function of the mass. The ionizing flux may then eventually 
overwhelm the accretion of neutral material, and this does 
play a role in our simulation. The third way in which a 
swamped source may still ionize some of its surroundings 
occurs if the accretion flow is anisotropic, so that radiation 
may escape from the source in some directions. It is this 
third way that is of most importance here; the dense 
gas around our ionizing sources is highly anisotropic. In 
Figures 7 and 8, we plot at two different epochs (at the 
time ionization was switched on, and 1.62 Myr later). 
Hammer projections depicting, from the point of view of 
randomly-chosen ionization sources, the directions in which 
ionizing radiation is able to penetrate the surrounding gas 



to a depth further than 5pc. Areas coloured black are those 
in which the radiation is absorbed before reaching this 
radius. For all the sources, large fractions of the sky are 
opaque to their radiation. At the earlier epoch, the average 
fraction of sky opaque beyond 5 pc as seen by all sources 
is ^ 70%. At the later epoch, this fraction has dropped to 
~ 45%, so clearly the effect of ionization increases as time 
passes, but rather slowly. 

Since the accretion flows are too vigorous to be 
disrupted by the ionizing radiation and they cover a large 
fraction of the sky as seen form most of the ionizing sources, 
the effect of feedback on most of the gas in the cluster, and 
hence on the star formation process is small. We terminated 
the simulation after ionization has been active for ~ 2.2Myr, 
at which point the age of the system is ~ 6Myr and the 
first supernova to occur is due in ~ O.SMyr, at which point 
the first 0-stars will reach 3Myr in age and leave the main 
sequence. By this time, 2 x 10* M© of gas was unbound, 
about two of thirds of it ionized, representing only ^ 2% of 
the system's total mass. This corresponds to a rate of mass 
loss of ~ lO~^M0yr~^. Given that the cluster contains on 
average 10s to a few hundred 0-stars, this rate of mass-loss 
is rath er low compared to those reported by, e.g.. iMazurekl 
(|l980l ), due to the quenching effect on the feedback of the 
accretion flows onto the ionizing sources. Although we have 
only allowed ionization to act for ^ 2.2 Myr, we think it 
highly unlikely that the effects of ionization will change 
significantly in the remaining time until supernovae begin 
to explode in the protocluster. 



4 DISCUSSION AND CONCLUSIONS 

We find that in large stellar systems whose formation has 
been triggered by the decay of a turbulent velocity field, 
star formation produces a complex hierarchical structure 
with clustering on many different length and mass scales. 
Star formation generally follows the dense filaments gener- 
ated by the velocity field, with the most massive and densest 
clusters found at the junctions of the filaments. Conversely, 
the protocluster exhibits large voids in which no star forma- 
tion occurs at all. The freefall time of our model system is 
~ 4Myr, so that dynamical interactions are unable to erase 
the large-scale structure of the protocluster on the timescale 
on which stellar feedback is expected to act and in partic- 
ular cannot greatly influence such structure before the first 
supernovae begin to detonate. The configuration of stellar 
clusters and gas produced by the turbulent velocity field 
therefore defines the backdrop against which feedback oper- 
ates. 

Photoionizing radiation has little impact on the evo- 
lution of our model protocluster and in particular is unable 
to disrupt the dense gas in which most of the star formation 
is occurring due to the presence of strong accretion flows 
onto the ionizing sources. As a result, photoionizing feed- 
back on this system is restricted to filling pre-existing voids 
in the cold gas with hot ionized material and the principal 
dynamical effect of feedback derives from the pressure ex- 
erted by this hot gas on the neutral gas. We compute at 
several times the average density of the ionized gas, which 
would have a temperature of 10'* K and select all the neutral 
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(a) Column density plots of cold neutral gas and stars (left panel), cold neutral gas and hot ionized gas (centre panel) and hot ionized 
gas alone (left panel) 0.66 Myr after ionization was turned on. 
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(b) Column density plots of cold neutral gas and stars (left panel), cold neutral gas and hot ionized gas (centre panel) and hot ionized 
gas alone (left panel) 1.08 Myr after ionization was turned on. 
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(c) Column density plots of cold neutral gas and stars (left panel), cold neutral gas and hot ionized gas (centre panel) and hot ionized 
gas alone (left panel) 2.18 Myr after ionization was turned on. 



Figure 6. Evolution of the ionized and neutral gas. 
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Figure 7. Hammer projections showing directions in which ionizing radiation is absorbed before reaching a radius of 5pc (black areas) 
from the point of view of three sources at the time when ionizing radiation was switched on. 
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Figure 8. Hammer projections showing directions in which ionizing radiation is absorbed before reaching a radius of 5pc (black areas) 
from the point of view of three sources 1.62 Myr after ionizing radiation was switched on. 



10° 



10- 



10- 



10- 



-1 1 1 1 r- 



-I 1 r- 




iD~' uf id' m* 



if 



Figure 9. Comparison of the initial (j;— axis) and final (x— axis) 
densities of randomly— selected neutral (black) and ionized (red) 
particles. 
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Figure 10. Plot of cumulative mass against initial density for 
the neutral (black) and ionized (red) gas. 



particles whose thermal pressure is in excess of the average 
thermal pressure of the ionized gas. The particles so selected 
are then too dense to be affected by the pressure of ionized 
gas. We find in this way that only a small amount, ~ 3% of 
the gas in the cluster at any time is too dense to be affected 
by the low-density ionized gas produced by the 0-stars. 
We illustrate this qualitatively in Figure [TT] where we plot 
the cumulative mass against density for neutral gas in black 
and ionize gas in red 1.62Myr after ionization was turned on. 
The vertical dashed line is the median density of the ionized 
gas and the vertical dotted line is the densest neutral gas 



whose pressure is less than that in ionized gas at the me- 
dian density, assuming all neutral gas has a temperature of 
40K. Neutral gas to the right of the dotted line would then 
be too dense to be strongly affected by the thermal pressure 
in the ionized gas. 

In principle then, the ionization should be able to signif- 
icantly influence the evolution of the protocluster's gas re- 
serves. In Figure 1121 we plot the positions of the gas whose 
density is too large to be affected by the HII. We neglect 
here the ram-pressure of the cold gas against which the HII 
must also contend, since the gas is infalling towards the ion- 
izing sources. The dense gas clearly (and not surprisingly) 
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Figure 11. Cumulative mass plotted against density for 
randomly— selected neutral (black) and ionized (red) particles 1.62 
Myr after ionization was turned on. The dashed vertical line in- 
dicates the median density of the ionized gas and the vertical 
dotted line indicates the densest cold gas whose pressure would 
be less than that of the median-density ionized gas if the cold gas 
were all at 40K. Neutral gas to the right of the dotted line cannot 
therefore be strongly influenced by pressure from ionized gas. 



correlates strongly with the positions of the clusters in the 
simulation - the stars are to be found embedded in the dens- 
est gas. In order to significantly affect the evolution of the 
cluster, photoionization would need to clear this material 
away. Since it can do so only very slowly, as shown in Fig- 
ures [7] and [S] the hot gas instead leaks into the low-density 

voids. 

Classical mode ls (e.g lBodenheimer et al.|[l979l : iLarsorJ 

1 19871 : lMatzne3l2002l ) indicate that ionization-driven cham- 
pagne flows should be an efficient mechanism for disrupting 
molecular clouds and t his is very likely to be true in many 
systems, such as Orion (|Hennev et aLlfioOSl ). where 0-stars 
are located near the edges of clouds. However, in cases where 
O-stars are embedded in highly inhomogeneous and often 
filamentary gas deep in cloud interiors, photoionization may 
not be so e ffective. In this work and that of several other 
authors fe.g lYorke et~aLlll989l : [Dale et akllioosl : iPeters etld] 
l2010al lbh. embedded O-stars struggle to disrupt their natal 
clusters by photoionization. This does not imply that this 
must always be the case, but it does suggest that gas expul- 
sion by ionizing radiation is more difficult than implied by 
earlier models. This difference emerges principally because 
strong accretion flows onto the ionizing sources suppress the 
growth of HII regions, a process which was previously often 
neglected. 

; Comeroi] 



S tellar winds may assist the UV radiation (cl , 

1 19971 : ICapriotti fc K ozminski' '200l'; 'Freve r et al.l l2003l l. al- 
though we found in , Dale fc Bonnell (2008) that cluster dis- 
ruption by winds alone may be a rather slow process. The 
oldest O-stars in our simulation will soon explode as su- 
pernovae and it is possible that the mechanical and ener- 
getic feedback from the few hundred O-stars that reside in 
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Figure 12. Positions of the Ri 3% of the gas too dense to be 
influenced by the thermal pressure of ionized gas at the time 
when ionization is turned on. 



our clusters would succeed where photoionization has failed. 
iTenorio-Tagle et al.l l|l985l ). for example, concluded that a 
single supernova could disrupt ~ 10* of molecular cloud ma- 
terial. The gravitational binding energy of our model system 
is a few X 10^^ erg and it has formed the equivalent of a few 
hundred O-stars by an age of 6Myr, so that supernovae 
should in principle have no trouble un binding it. Alterna - 
tively, as suggested and investigated bv lClark et al.l l|2005h , 
molecular clouds and their embedded clusters (or cluster 
complexes) may not be bound in the first place. 
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